#include<iostream>
using namespace std;
#include<cmath>
#include<fstream>

#include"library.h"

void derivative(int n, double x, double *y, double *dydx){
	dydx[0] = y[1];// dx/dt = v
	dydx[1] = - y[0];// - 0.1 * y[1];// + 0 * sin(1.1*x);// dv/dt = -x - gamma v + F_0 sin(omega t)
}

int main(){
	int i,j, n=2, N=10000;
	double a=0, b=10;
	double h = (b-a)/N;
	cout<<"h="<<h<<endl;
	double *y = new double [n];
	double *ynext = new double [n];
	ofstream fp("classicHO.txt");
	y[0]=1; y[1] = 0;
	fp<<a<<"\t"<<y[0]<<"\t"<<y[1]<<"\n";
	for(i=0;i<N;i++){
		rk4(a+i*h, h, n, y, ynext, derivative);
		for(j=0;j<n;j++) y[j] = ynext[j];
		fp<<a+(i+1)*h<<"\t"<<ynext[0]<<"\t"<<ynext[1]<<endl;
	}
	fp.close();
	return 0;
}
